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^ Abstract 

The objective of this article is to apply recent developments in geo- 
, metric optimal control to analyze the time minimum control problem of 

' dissipative two-level quantum systems whose dynamics is governed by the 

Lindblad equation. We focus our analysis on the case where the extremal 
' _i ' I Hamiltonian is integrable. 
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1 Introduction 

' We consider a dissipative two-level quantum system whose dynamics is governed 

CO I by the Lindblad equation which takes the following form in suitable coordinates 
q — {x, y, z), i.e., in the coherence vector formulation of density matrix [U I16j : 

'NT ' 

0> ■ X — —Tx + U2Z 

^ ■ y = -Ty~uiz (1) 

■ z — J- — 7+z 4- uiy — U2X. 

We refer to (TU] and [T7] for the details of the model. We recall that x and y 
' are related to off-diagonal terms of the density matrix of the system and z to 

the difference of population between the two states. In ([T]), A — (r,7_,7+) is 
a set of parameters such that T > 7-|_/2 > and 7+ > |7_| which describes 
the interaction of the two-level system with the environment. More precisely, 
r is the dephasing rate and 7+ and 7_ are respectively equal to 712 + 721 and 
712 — 721 where the coefficients 712 and 721 are the population relaxation rates. 
The control is the complex Rabi frequency u = ui + iu2 of the laser field which 
is assumed to be in resonance with the frequency of the two- level system 
The physical state belongs to the Block ball, l^l < 1, which is invariant for the 
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dynamics considered. If they are many articles devoted to optimal control of 
quantum systems in the conservative case (see e.g. jl2^), the dissipative case is 
still an open problem. 

The system can be written shortly as a bilinear system 

q ^ Foiq) + UiFiiq) + U2F2iq) (2) 

and in order to minimize the effect of dissipation, we consider the time minimum 
control problem for which, up to a rescaling on the set of parameters A, the 
control bound is |u| < 1. The energy minimization problem with the cost 

\u\'^dt where the time T is fixed but the control bound is relaxed can also be 
considered and shares similar properties. 

A first step in the analysis of such systems is contained in [17) . Assum- 
ing u real, the problem can be reduced to the time-optimal control of a two- 
dimensional system 

y = ~Ty-uiz (3) 

i = 7_ - 7+z -I- uiy, 

with the constraint < 1. For such a problem, the geometric optimal control 
techniques for single-input two-dimensional systems presented in [13 succeed to 
make the time-optimal synthesis for every values of parameters (r,7„,7+). 

In order to complete the analysis in the bi-input case, a different methodology 
has to be applied and we shall make an intensive use of techniques and results 
developed in a parallel research project to minimize the transfer of a satellite 
between two elliptic orbits (see [lO]). Such techniques are two-fold. 

First of all, the maximum principle will select extremal trajectories, can- 
didates as minimizers and solutions of an Hamiltonian equation. A geometric 
analysis will identify the symmetry group of the system and find suitable co- 
ordinates to represent the Hamiltonian. A consequence of this analysis is the 
fact that, in the case 7_ = 0, the extremal system is integrable and if F = 7+, 
the problem can be in addition reduced to a 2D- almost Riemannian problem 
on a two-sphere of revolution for which a complete analysis comes from [B] . We 
take advantage of this property to analyze the general integrable case 7_ = 
using continuation methods on the set of parameters, while the analysis fits 
in the geometrical framework of Zermelo navigation problems [3]. This case 
corresponds to the physical situation where the population relaxation rates are 
equal. The mathematical tools presented in this article and a few numerical 
simulations are sufficient to complete the analysis for 7_ = 0. The generic case 
where 7_ ^ is treated in a forthcoming article [5] combining mathematical 
analysis and intensive numerical computations. 

Secondly, having selected extremal trajectories, second-order conditions us- 
ing the variational equation and implemented in the cotcot code allow to 
determine first conjugate points forming the conjugate locus which are points 
where extremals cease to be locally optimal. Combined with the geometric anal- 
ysis, we can construct the cut locus which is formed by points where extremals 
cease to be globally optimal. 

The organization of this article is the following. In section [21 we recall the 
maximum principle and the concept of conjugate points associated to second- 
order optimality conditions. In section [31 we present the geometric analysis 
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of the system, followed in section 0] by a thorough analysis of the so-called 
Grusin problem on a two-sphere of revolution. This problem is generalized into 
a Zermelo navigation problem suitable to our analysis. In section we study 
the properties of the extremals to analyze the optimal trajectories, combining 
analytic and numerical methods. 



2 Geometric optimal control 
2.1 Maximum principle 

We consider the time minimum problem with fixed extremities qo and qi . For a 
smooth system written q = F{q, u), with q S M" and a control domain U which 
is a compact subset of M'", we have: 

Proposition 1. // {q,u) is an optimal control trajectory pair on [0,T] then 
there exists an absolutely continuous non-zero vector function (j),po) G M" x M 
such that almost everywhere on [0,r], we have: 

■ 9H. , . dH. . 

and 

Hiq,p,u) ^ Miq,p) (5) 

where H{q,p,u) — {p, F{q,u)) + po, Po being a non positive constant; M is 
defined by M{q,p) = MaXueuH{<lTPiU) and M is zero everywhere. 

Definition 1. The mapping H is called the pseudo-Hamiltonian. A triple 
{q,p,u) solution of ^ and is called extremal while the component q is 
called an extremal trajectory and p is an adjoint vector. 

Application: We consider the time-minimum control problem for a system 
of the form q = Fo{q) + '^iFi{q) with m > 2, u — (ui, ■ ■ ■ , Um), \u\ < 1. 

We introduce the Hamiltonian lifts Hi — {p, Fi{q)), i = 0,1, ■ ■ ■ , m of the vector 
fields Fi and the set S such that Hi — for i = 1, ■ ■ ■ ,m. Then the maximization 
condition ^ leads to the following result. 



Proposition 2. Outside S, an extremal control is given by Ui ~ Fli / y'Y^'^^ Hf , 
i = ,m and extremal pairs z — {q,p) are solutions of the smooth true 

Hamiltonian vector field H r{z) with Hr{z) — Hq{z) + (X^I^i -^fi^))^^'^- 

Definition 2. The surface S is called the switching surface and the solutions 
of Hr{z) are called extremals of order zero. To be optimal, they have to satisfy 
Hr{z) > and those with Hr{z) = are called abnormal. 

An important but straightforward result is the following proposition [8] . 

Proposition 3. Extremal trajectories of order zero correspond to singularities 
of the end-point mapping 

E'i^''^ ■.ueL°°[0,T]^^q{T,qo,u) (6) 



where q{-, qo, u) denotes the response to u with initial condition qo and such that 
the control is restricted to the (to — l)-sphere |u| = 1. 
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2.2 Second-order optimality conditions 

From proposition [21 we can apply the concepts and algorithms presented in [7] 
to compute second-order optimality conditions in the smooth case, when the 
control domain U is a manifold, u being restricted to the unit sphere. The 
framework of this computation is next recalled. 
The concept of conjugate point: 

Since J7 is a manifold, we may assume locally that U = M™^^ and the 
maximization condition ^ leads to 

du ' du^ ~ 

Our first assumption is the strong Legendre-Clebsch condition: 
(HI) The Hessian d^H/dv?' is negative definite along the reference extremal. 
From the implicit function theorem, an extremal control can be locally defined 
as a smooth function of z = {q^p) and plugging u into H defines a smooth true 
Hamiltonian Hr- 

Setting Af — R" and using Hamiltonian formalism, we introduce: 

Definition 3. Let z — {q,p) he a reference extremal defined on [0, T]. The 
variational equation 

5z = dHr{z{t))5z (8) 

is called the Jacobi equation. A Jacobi field is a non trivial solution 5z ~ 
{5q,5p). It is said to be vertical at time t if Sq{t) — 0. 

The following standard result is crucial. 

Proposition 4. Let Lq be the fiber T*^M and let Lt = exp[tiJr(io)] be its 

image by the one parameter subgroup generated by H r- Then Lt is a Lagrangian 
manifold whose tangent space at z{t) is spanned by Jacobi fields vertical att = 0. 
Moreover, the rank of the restriction to Lt of the projection H : (q,p) i— > g is at 
most [n — 1). 

We next formulate the relevant generic assumptions using the end-point map- 
ping. 

Assumptions: 

1. (H2) On each subinterval [to,^i], ^ < to < ti < T^ the singularity of 
^9(to),ti-to jg Qf codimension one for it|[tQ j^j. 

2. (H3) We are in the normal case Hr ^ 0. 

As a result, on each subinterval [to,^i] there exists up to a positive scalar an 
unique adjoint vector p such that {q^p, u) is extremal. 

Definition 4. We fix qo — q{0) and we define the exponential mapping; 

expg^ : (p(0),t) ^ H(exp[ti?,((?(0),p(0))]) 
where p{0) is a {n — 1) dimensional vector, normalized with iJ,. = f . 
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Definition 5. Let z = {q,p) be the reference extremal on [0,T]. Under our 
assumptions, a time < tc < T is called conjugate if the mapping exp^^^ is not 
an immersion at (p(0), tc) and the point q{tc) is said to be conjugate to go- M^e 
denote tic the first conjugate time and C{qa) the conjugate locus formed by the 
set of first conjugate points considering all extremal curves. 

We get: 

Theorem 1. Let z{t) = {q{t),p{t)) be a reference extremal on [0, T] satisfying 
assumptions (HI), (H2) and (H3). Then the extremal is optimal in the L°° - 
norm topology on the set of controls up to the first conjugate time tic- Moreover, 
if t ^ q{t) is one-to-one then it can be embedded into a set W , image by the 
exponential mapping exp^^Q^ of N x [0, T], where N is a conical neighborhood of 
p(0). ForT < tic, the reference extremal trajectory is time minimal with respect 
to all trajectories contained in W . 

In order to get global optimality results, it is necessary to glue together such 
micro-local sets. We need to introduce the following concepts. 

Definition 6. Given an extremal trajectory, the first point where it ceases to 
be optimal is called the cut point and taking all the extremals starting from q^, 
they will form the cut locus Cut{qo). The separating line SL{qQ) is formed by 
the set of points where two minimizers initiating from go intersect. 

3 Geometric analysis of Lindblad equation 

3.1 Symmetry of revolution 

If we apply to system ([T]) a change of coordinates defined by a rotation of angle 
6 around the z— axis: 

X ~ X cos 9 -\- ysin9 
Y = —xsm9 + ycos0 
Z = z, 

and a similar feedback transformation on the control: 

vi = ui cos 9 U2 sin 9, V2 — ~ui sin + 1*2 cos 9, 
we obtain the system 

X = -TX + V2Z 

Y = -TY-viZ (9) 
Z = "/^ — ^^Z + viY — V2X. 

Hence, this defines a one dimensional symmetry group and by construction 
|u| = \v\. Therefore, we deduce that the time-minimum control problem (or the 
energy minimization problem) are invariant for such an action. Using cylindric 
coordinates 

x — rcos9, y — rsm9, z^z 
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and the dual variables p — (pr,p9,Pz), the Hamiltonian Hr takes the form 

Hr = Ho + iHl + Hlf^ 

= {-Vrpr + (7- - 7+)Pz) + {z'^vl + + 'r^v\ - 4,zrprPz)^^^ . 

In particular, the Bloch ball is foliated by meridian planes 6 — constant in 
which the time-minimum synthesis is the one associated to system ([3]), where 
the control is scalar and described in [17] . More precisely, we have: 

Proposition 5. For the time-minimum control, 6 is a cyclic coordinate and pg 
is a first integral of the motion. The sign of 9 is given by pg and if pg = then 
9 is constant and the extremal synthesis for an initial point on the z-axis is up 
to a rotation given by the synthesis in the plane 9 = 0. Up to a rotation, the 
control u can also be restricted to the single-input control {ui,Q). 

Proof. The proof is a generalization of the geometric situation encountered in 
[5]. For the Hamiltonian vector field Hr, the points on the z-axis correspond 
to a polar singularity and the extremals starting from the z-axis are contained 
in meridian planes 9 = 9{Q). Hence, pg is constant and extremal curves in the 
plane 9 = 9{Q) are solutions of system ([3]). □ 

3.2 Spherical coordinates 

More properties can be seen using spherical coordinates: 

X = p sin (j) cos 9, y = p sin sin 9, z = p cos 

and a similar feedback transformation. We obtain the system: 

p = 7_ cos (j) — p{j+ cos^ (f> + T sin^ (j)) 
■ 7- sin , sin(20) 

(j) = + — 7+-r) + i;2 (10) 

P 2 

9 = — (cot0)wi. 

and the corresponding Hamiltonian 

Hr = [7- cos (j) — /o(7+ cos^ + r sin^ 4>)]Pp 

, 7_sin0 sin(2(7!)) , rz. Z. ^ 

+[-^^ + —^(7+ - r)W + \jv\ +pI cot" 0. 

In this representation, {(j), 9) are the spherical coordinates on the unit sphere of 
revolution around the z-axis: 9 is the angle of revolution and (j) G]0,7r[ is the 
angle of the meridian, </) = 0, tt correspond respectively to the north and south 
poles. 

3.3 Lie brackets computations 

In order to complete the analysis, we immediately compute the Lie brackets up 
to length 3 for the system written in Cartesian coordinates as: 

q = (Gaq + vq) uidq + 
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where the Gi's are the matrices 

/ -r \ 

Go = -r , Gi = I -1 1 , G2 



-r 











-r 











-7+ 





and *wo is the vector (0,0, 7_). It can be lifted into a right-invariant control 
system on the semi-direct product GL(3,Ili) xglR^ identified to the subgroup of 
matrices of GL{4, M) of the form: 



1 
V g 



, v€W^,ge GL{3, 




acting on the subset of vectors of : ^ ^ J , g G M'^ To construct affine vector 

fields, we use the induced action of the Lie algebra {a, A) ■ q = Aq + a and Lie 
brackets are given by 

[(a, A), {b, B)] = [Ab - Ba, AB - BA]. 

The control distribution is £> = Span{Gi, G2} and we have: 



[Gi, G2] = G3 



In particular, we obtain that {Gi,G2}a.l. = so{3) and hence the system on 
SO{3): 

^ = {uiGi+U2G2)X 

is controllable. For the linear action, it defines a controllable system on the unit 
sphere. This action has however singularities: 

• at 0, the orbit is 0. 

• the set on where Gi and G2 are collinear is the whole plane z = and 
restricted to the unit sphere of revolution, it corresponds to the equator. 

To analyze the effect of the drift term associated to dissipation, we use 



[Go,Gi] = (r-7+) I 1 I , [Go,G2] = (7+-r) 
1 

and [Go, G3] = 0. Moreover, we have 

\ / 1 

[Gi,[Go,Gi]] =2(r-7+) I 1 I , [G2,[Go,G2]] = 2(7+-r) I 

0-1 

Those computations reveal the singularity at 7-)- = F, 7_ = that we describe 
in the next proposition. 
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Proposition 6. In the case 7_ = 0, 7+ = T , the radial component p is not 
controllable and the time-minimum control problem is an almost Riemannian 
problem on the two-sphere of revolution for the metric in spherical coordinates 
g — dcjp' + tan^ 0d6'^ with Hamiltonian H = ^ (p| + pg cot^ <j)) . 

Proof. If 7_ and 7+ = F then the first equation of (fTU)) becomes p ~ — Fp 
and p is not controllable. Hence, the time-minimal control problem reduces to 
the problem of controlling (0, 9) in minimum time where the associated true 

Hamiltonian is y^P^ + Pq cot^ 0. The time minimum problem is equivalent to 

minimizing the length for the metric d(f>^ + tan^ (pdO^ . According to Maupertuis 
principle, we can replace the length by the energy with corresponding Hamilto- 
nian i (p^ -f pg cot^ (j)) . □ 

Definition 7. The almost Riemannian metric g = d(j)^ -|-tan^ (j)dO^ is called the 
Grusin model on the two-sphere of revolution. 

Such a metric appears in quantum control in the conservative case jl2j and a 
similar metric is associated to orbital transfer [4]. It will be analyzed in details 
in section [4] since it is the starting point of the analysis in the general case using 
a continuation method on the set of parameters. 

Another consequence of the previous computations is the controllability 
properties of the system and the structure of extremal trajectories. 

3.4 Controllability properties 

We recall that the Bloch ball \q\ < 1 is invariant. Indeed, introducing p^ — jqp, 
we get: 

pp^~T{x^ +y'^)--i+z'^ + -t^z <Q (11) 

which is strictly negative on the unit sphere except ii -\- = Q, \z\ = \ and 
7+ — |7-|. Using the representation PT7| of the system in spherical coordinates, 
it is clear that we can control the angular variables 4> and 9 if the controls are 
not uniformly bounded. If |u| < 1 then we have restrictions depending upon 
the set of parameters. 

For 7_ = 0, the system is homogeneous and g = is a fixed point. The 
accessibility set in fixed time is with non empty interior for a non-zero initial 
point, except in the case 7-1- = F which corresponds to the Grusin model and 
for which the time and energy minimization problem are equivalent. 
The controllability properties for |m| < 1 are clear in this case. Indeed, if 
I7+ — F| < 2 then we can compensate the drift by feedback for the system on 
the two-sphere of revolution, while it is not the case for |7_|_ — F| > 2. 

Proposition 7. Let qq and qi be two points in the Bloch ball \q\ < 1 such that 
qi is accessible to go- Then there exists a time-minimum trajectory joining qo 
to qi . Moreover, every optimal trajectory is 

1. either an extremal trajectory with pe = 0, contained in a meridian plane, 
time-optimal solution of the two-dimensional system (0] where u = (ui, 0). 

2. either connection of smooth extremal arcs of order 0, solutions of the 
Hamiltonian vector field H r with pg ^ 0, while the only possible con- 
nections are located in the equatorial plane <j) — 7r/2. 
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Proof. The control domain is convex and the Bloch ball is compact. Hence, we 
can apply the Filippov existence theorem |15| . In order to get a regularity result 
about optimal trajectories, much more work has to be done. This is due to the 
existence of a switching surface S : Hi = H2 = in which we can connect 
two extremals arcs of order 0, provided we respect the Erdmann-Weierstrass 
conditions at the junction, i.e., the adjoint vector remains continuous and the 
Hamiltonian is constant. The set E can also contain singular arcs for which 
Hi = H2 = holds identically. Hence, we can have intricate behaviors for such 
systems. In our case, the situation is simplified by the symmetry of revolution. 

Indeed, ii pe — then the singularities are related to the classification of 
extremals in the single-input case, which is described in [17j . We cannot connect 
an extremal with pg ^ where pg is the global first integral xpy — ypx to an 
extremal where pg ~ Q since the adjoint vector has to be continuous. 

Hence, the only remaining possibility is to connect two extremals of order 
with pg ^ Q aX a, point of S leading to the conditions p^ = and pg cot (/) — in 
spherical coordinates. Since pg ^ 0, one gets (p — tt/2. The result is proved. □ 

Remark: The classification of extremal trajectories near the equatorial 
plane is described in proposition 1191 

4 The Grusin model on a two-sphere of revolu- 
tion with generalizations to Zermelo naviga- 
tion problem 

The Grusin model g — dcfP' + tan^ cjidO'^ is a special case of metrics of the form 
dcj)^ + G{(t))d9'^ on a two-sphere of revolution such that: 

• (HI) G'(0) 7^0on]0,7r/2[. 

• (H2) G{ti ~ (j)) — G{(p) (reflective symmetry with respect to the equator). 

They appear in optimal control in the orbit transfer, smooth at the equator or 
with a polar singularity, in quantum control and in various geometric problems, 
e.g., Riemannian problems on an ellipsoid of revolution. The importance of 
this control problem has justified the recent analysis of [Ij that we complete 
next, using Hamiltonian formalism in order to make generalizations. We first 
interpret the Grusin model as a deformation of the round sphere. 

Definition 8. The standard homotopy between the Grusin model and the round 
sphere on the two-sphere of revolution is g\ — d(p'^ + G\{X)d6^ where G\{X) = 

i^^X ' = ^^^^ ^ "-^^ ^ t*^' -"^l- 

By construction, the metric is analytic for A G [0, 1[ and for A = 1, we have 
the Grusin model with a pole of order 1 at the equator. 

The first objective of this section is to show stability results concerning such 
metrics. We have the following general result [S]. 

Proposition 8. Let d<p'^ + G{(l))d0^ be a smooth metric on a two-surface of 
revolution. Then, 

1. Extremals are solutions of the Hamiltonian H = ^{p'^ + c^^j) '^'^'^ ^'^c- 
length parametrization amounts to restrict to H = 1/2. 
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2. If ip is the angle of an unit-speed extremal with a parallel then pe = 
\fG cos tp is a constant and the extremal flow is Liouville integrable with 
two commuting first integrals H and pg . 

3. The Gauss curvature is K = — ^ ■ 

We next make a complete analysis of the family of metrics g\. 



4.1 Curvature analysis 

We have the following proposition. 

Proposition 9. For the family of metrics g\, we have: 
1. The Gauss curvature is Kx = ''^^li^xsin^^^)'^ ■ 
^- K'{<t>) = (i^rf ^4 [5(1 - A) - 4Acos^ <!>]. 

Hence K{(j)) is non-constant and monotone non decreasing from the north pole 
to the equator for A g]G, 1/5], while for A g]1/5; 1[ it admits a minimum. For 
A e]0, 1[, the curvature is maximum on the equator. The hmit case A = 1 
corresponds to the Grusin case, for which the curvature is negative everywhere 
and tends to — oo when cj) tends to 7r/2. 



4.2 Geometric properties 

We next present the main properties of the extremal flow for a metric on a 
two-sphere of revolution g = dcjp' + G{(j))d9^ satisfying (HI) and (H2) where g is 
smooth, except may be at the equator where it can admit a pole of order one. 

2 

For such a family, we consider the smooth Hamiltonian H = i(p^ + 
and we restrict extremal curves to the level set H = 1/2. Fixing pe, the param- 
eterized family of corresponding Hamiltonians described the evolution of the 

2 

(p0, (j)) variables as solutions of a mechanical system for which V{4i) = 
plays the role of potential. For pg = 0, we get the meridian solutions. Hence, 
we can assume pg ^ Q. Using assumption (HI), the only equilibrium point is 
for 4> = 7r/2 and p^ = 0. This leads to the equator solution in the regular case. 

For the remaining trajectories, the level set iJ = 1/2 is sufficient to analyze 
the behaviors of </). Indeed, it is a compact set, symmetric for the two reflections 
with respect to the 0-axis and the equator (/> = 7r/2 and defined respectively by 
the two transformations: i— > -~P4, and (j) i-^ tt — <f). Every trajectory is periodic 
and 4' = 7r/2— ^ oscillates periodically between ipmax and —ipmax- There is also a 
relation between the period of oscillation T and the amplitude tjjmax, depending 
upon pg. 

By symmetry, every trajectory is defined by its restriction to a quarter of 
period, that is the sub-arc starting from the equator tp = and reaching tpmax- 
The trajectory starting from {(j){0) , p^{Oj) and reaching tt — (j){0) after passing 
i^max corresponds to a point rotating on the level set H = 1/2 and is chased by 
a point associated to the trajectory starting from (^(0), — ^^(0)) and reaching 
TT — 0(0). They are distinct if ^^(0) ^ 0. Moreover, using the assumption 
(H2), we deduce easily that for fixed pg, the extremals starting from {4>{0), 9{0)) 
with respectively p^(0) and — p^(0) intersect with equal length on the antipodal 
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parallel. They are distinct if p,^(0) 7^ 0. The case ^^(O) = corresponds for 
an initial condition not on the equator to tangential arrival and departure at 
parallels 0(0) and tt — 0(0); ^^(0) = gives the equator solution in the non 
singular case. For more details see [S]. 

The only difference in the singular case is that the equator is not solution, and 
for trajectories departing from the equator the extremals are always tangential 
to the meridian, while the first return to the equator can be arbitrarily closed 
from the initial point. 

Finally, another obvious symmetry is a reflectional symmetry with respect 
to the meridian obtained by changing pe into —pg. 

As a consequence of this analysis, we deduce: 

Proposition 10. Let g = dcfy^ + G{4>)d6'^ be a metric on a two-sphere of revo- 
lution, satisfying (HI) and (H2) and smooth except may he at the equator where 
it can admit a pole of order one. 

1. Then except the meridian and the equator solution in the regular case, 
every extremal is such that ip = 7r/2 — oscillates periodically between two 
symmetric parallels. The first return mapping to the equator is 

R:pe e]0, v/G(V2)H A%e), 
where AO is the corresponding 9 — variation of the extremal. 

2. Assume pg ^ and Pc/){0) ^ then 

a) Fixing pg, changing p^(0) into —p^{0) gives two distinct extremals with 
equal length intersecting on the antipodal parallel. 

h) Fixing p^(0) and changing pg into —pg gives two distinct extremals with 
equal length intersecting on the opposite meridian. 



4.3 Integrability 

For the family of metrics gx which fit in the previous geometric framework, we 
can be more precise and make a complete analysis. For a fixed value of A, the 
Hamiltonian is: 

and corresponds for A = 1 to the Grusin case. Using = p^, we get: 

lr;2 , ^sin^<^)i 



which can be written 



Hx~ > .2 

^ sm ( 



i^A = ^[0'+P^(cot20 + l-A)]. 



Therefore, Hx = Hi + ^pj{l — A) and parameterized by arc-length: Hx = 1/2, 
one gets the level set Hi — ^ — ^Pg{l — A). 

Hence the integration of the Grusin case gives the general solution, and from 
the homotopy, the corresponding extremals fit not only in the same geometric 
framework, but also have the same transcendence. 
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Lemma 1. The family of Hamiltonians H\ admits two first integrals in invo- 
lution for the Poisson bracket ( independent of X) pg and Hi = ^ [p^ + pg cot^ (j>] . 

We next outline the integration method in the Grusin case to provide the 
computation of the first return mapping to the equator R obtained in [5]. We 
have Hi — ^{(f^ + vcoi^ (f>), v > and fixing the level set to 1/2, we get: 



(^) 



2 



1 — {v + 1) COS 



2 



dt sm" 

Taking the positive branch, we must evaluate the following expression 



=t. 



(1 - 1)C0S2 </))l/2 

To integrate, we use the relation 

COS (bd(h 1 



■ arcsm m sm ( 



\/l - m2 sin^ (j) m 
to deduce the form of the component (f> of the general solution: 

(f>(t) = arcsinf— (sin(mt + K))] + tt/2. 
m 

To complete the integration, we write: 

A Pe 



■ 2 J. ^PS 

sm 6 



and we use the formula: 



f dx 1 n , 

/ 7^ — = arctan[v 1 — atana;J 

J 1 — a sin X \J\ — a 

for a < 1 with the relation cos (\) — ~ — sin(mt + K\ A straightforward compu- 
tation then leads to B(t). 

Proposition 11. For the family of metrics g\, we have: 

anpe X 



In particular, if a > then R'{pe) < < R"{pe) on ]0, ^JG{t:/2)[. 

This property allows to evaluate conjugate and cut loci for the family of 
metrics that we next describe [6]. 

4.4 Conjugate and cut loci 

We have: 



Theorem 2. • For A = (round sphere), the conjugate and cut loci of any 
point are reduced to the antipodal point. 
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• For < A < 1, the conjugate locus of a point different from a pole is 
diffeomorphic to a standard astroid, while the cut locus is a single branch 
of the antipodal parallel. Both are symmetric with respect to the opposite 
meridian. 

• For A = 1 (Grusin case), the conjugate and cut loci of a point different 
from a pole and not on the equator are as above. For a point on the equator, 
the cut locus is the equator minus this point and for the conjugate locus, 
the cusps on the equator are transformed into folds at this point minus 
this point. 

Geometric interpretation: 

For the class of metrics g\, the situation is clear. For the round sphere, all 
extremals starting from the equator intersect at the same antipodal point and 
the first return mapping is constant. For < A < 1, the first return mapping 
is monotone, and in the singular case R(j)e) as pe ^ +oo. Since the cut 
locus of a point of the equator is formed by intersections with the equator of 
symmetric extremals, in the homotopy, the cut locus is pinched into a point for 
A = 0, while it is stretched into the whole equator in the case A = 1. 



4.5 Zermelo navigation problem on the two-sphere of rev- 
olution 

We introduce the following definition for the Zermelo problem. 

Definition 9. A Zermelo navigation problem on the two-sphere of revolution 
is a time-minimum problem of the form: 

^==Fo(<z) + ^u,F,(g), \u\<l, 

2 = 1 

where the drift representing the current is of the form Fq (0) ^ + -Fq ('^) M 
Fi and F2 form outside the equator an orthonormal frame for a metric of the 
form g ~ dcj)^ + G{(f))d9^ . It is called reflectionaly symmetric with respect to the 
equator if 

• (HI) G'{(j)) on ]0,7r/2[ 

• (H2) G{tt -<P)^ G{(t>) 

. (H3) F^{^-cp) = -F^{<f), Fi = 0. 
It defines a Finsler geometric problem if \F()\ < 1 for the metric g. 
According to this classification, we have: 

Proposition 12. Assume 7_ = and consider the system IJOj) restricted to 
the two-sphere: 

. sin(20)(7+-r) , 

= h V2 

e = -{coicj))vi, \v\ < 1. 
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Then it defines a Zermelo navigation problem on the two-sphere of revolution 
where the current is Fq = ^i^^^^(7_|_ — T), the metric is g = d(fP + tan^ (pdO"^ 
with a singularity at the equator and the assumptions (HI), (H2) and (H3) are 
satisfied. The drift can be compensated by a feedback when I7+ — r| < 2, which 
defines a Finsler geometric problem on the sphere minus the equator. 

Controllability analysis 

The amplitude of the current is | sin(2^) (7_|_ — r)/2| and is maximum in the 
upper hemisphere for (p = 7i'/4, while it is minimum at the north pole and at 
the equator. Hence, more generally, we deduce the following proposition. 

Proposition 13. For I7+ — r| > 2, the current can be compensated in the north 
equator except in a hand centered at (f> ^ n /A, hence defining a Finsler geometric 
problem near the equator and near the north pole. 

The controllability analysis is straightforward and is related in the north 
hemisphere to the scalar equation: 

• sin(2V')(7+ - r) , 

^ ^ ^ ^'\^'+ L -V2, iP = tt/2 - 4), \V2\ < 1. 

Starting at ^ = with ti2 = — 1, to increase ip we meet a barrier corresponding 
to the singularity of the vector field. For instance, if 7+ — F > then we have 
a barrier when 1 = sin(2t/>)(7+ - F)/2, V e]0,7r/2[. 

5 The integrable case of two-level Lindblad equa- 
tions 

5.1 The program 

We now proceed to the analysis of the general case of a two-level Lindblad 
equation. The method is to start from the Grusin case and then to consider 
perturbations. This program succeeds only if 7_ = 0, leading to extremal flows 
described by a family of 2D-integrable Hamiltonian vector fields. 

5.2 The integrable case 

We observe that for 7_ = 0, the true Hamiltonian simplifies into: 

Hr = -p(7+ cos^ ^ + F sin^ ,p)pp + ^^^^^Mhl^lllp^ + y^P^Tpfc^t^^ 
and we immediately deduce: 

Proposition 14. For 7_ = 0, using the coordinate r = In p, the Hamiltonian 
takes the form: 

Hr = - (7+ COS^ (j} + T sir? (I})pr + ^^"(^'^)(^^+ ~ + y^p2 _^ p2 

Hence r and are cyclic coordinates andpr, pe are first integrals of the motion. 
The system is Liouville integrable. 
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We have the foUowing geometric interpretation. 

Proposition 15. For 7_ = 0, the Hamiltonian is associated if Pr < 
to the problem of minimization of r, while the case Pr > corresponds to the 
maximization of r, for |w| = 1. 

Proof. It is a consequence of the maximum principle. Another point of view is 
to consider the end-point mapping E : u i-^ q{t,qo,u). If u is restricted to the 
sphere |u| = 1 then the solutions of Hr parameterize the singularities of the end- 
point mapping. The case pr = corresponds to singularities of the end-point 
mapping for the system restricted to the two-sphere. In the extremum problem 
of r, with fixed time, Pr can be normalized to -1, or 1, while the level sets arc 
Hr = h. In the extremum problem of time, the Hamiltonian is normalized to 
or 1 for the minimum case, and or -1 for the maximum one. This is clearly 
equivalent by homogeneity. Hence, this gives a dual point of view. □ 

In order to indicate the complexity of the problem, we consider first the case 
of energy, which is equivalent to time from Maupertuis principle, in the Grusin 
case. 



5.2.1 The case of energy 

In the normal case, the true Hamiltonian is 

Hr = -Pr{l+ cos' + r sin^ 0) + Ei^^i^^ - D + i {pi cot' cP+pI). 

We fix the level set to h and using the relation p^ = (j) + EHfi^^p _ one 
gets \(j)' + V{(j)) = h where V{(j)) is the potential: 

y(^) = -pr{i+ cos' </. + rsin^ 0) - i^^^(r - 7+)^ + ^pI cot' 4>. 

To integrate, we use: 

d(t> 



±^2{h-V{4>)), 
and we must evaluate an integral of the form: 

dX 



where 

P(X)=4(l-X)[-X3(r-7+)2+x2(r-7+)(2p,+(r-7+))+X(2/i+2p,7++p2)-p2]^ 
with X = siv? (j). This corresponds to an elliptic integral. 



5.2.2 The time-minimum case 



In the time-minimum case, the computations of the extremal curves are more 
intricate because we cannot reduce the system to a second-order differential 
equation. The geometric framework is however neat because it is associated to 
a Zermelo navigation problem. 
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Computations: 



We set Q = yjP^ + Pg cot^ (j) a-nd the Hamiltonian is restricted to a level set e, 

where e = corresponds to the abnormal case and e = +1 to the normal case. 
This gives the following relation: 

-(7+ cos^ + rsin^ ^)p, + (7_ - p) Sm(^20) ^^ + Q = £ (12) 
and fixing pr and pg, the pair <p,p^ is solution of the system: 
^ = ^^sin(2^) + | 

z?^ COS d) 

P4> = (r-7+)sin(2(^K + (7+-r)cos(20)p^ + ^f-^. (13) 

Q sm (p 

Hence <j) — leads to: 

^^+~^^ sin(20)Q+p^ = O. (14) 
Using we deduce that is solution of a polynomial equation of degree 2: 

[(7+ - rf^-^^^ - 1] - (7+ - r) sin(20)[e + (7+ cos^ <j> + Tsin^ 0)p,,]p^ 
+ [£ + (7+ cos^ + r sin^ (/.)p^]2 - p2 cot2 = 0. (15) 
The discriminant of this polynomial is: 

A = 4[£ + (7+ cos^ (j) + rsin^ 0)Pr]^ +Pe cot^ ^[(7+ - r)^ sin2(20) - 4]. (16) 
From (fTi]) . we deduce: 

Plih. - r)^^ - 1] ^ -(7, - T^^pl cot^ ^. (17) 
Hence the set (0 = 0) n (i/r = e) is defined by the relation: 

[e + (7+ cos^ (t> + Tsin^ 0)pr]2 p2 ^^^2 _ *^'^+ ~ sin^ ,/)]. (18) 
Therefore, we have: 

Lemma 2. i. // (7+ - r)2Si£!lM - 1 =;z^ and A > then the level set 
Hr = s has two real roots p,f, which are distinct if A> 0. 

2. The intersection of (f> = with the level set Hr = s is given by A = 
which can be written: 

[e + (7+(i -x) + rx)prr = ^it^[i - (7+ - r^xii X)] (19) 

where X = sin^ (h. 
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Extremals analysis: 

From the previous analysis, we deduce that there are two types of extremal 
curves by considering the reduced system (jl3p describing the evolution of ((/>, p^). 
Compact case: 

It corresponds to the situation where the level sets Hr = e define compact 
surfaces in the 2-plane [(j),p^). In this case, if the reduced system is without 
singular point on the level set then the trajectory t ^ {(j){t),p^{t)) is a periodic 
trajectory with period T. 
The following lemma is clear. 

Lemma 3. The Hamiltonian Hr is invariant for the transformation {(f), p^) i—^ 
(tt- 0, -p^). 

A consequence of lemma |3] is the following. Assume that for (j)e,Pr) fixed, a 
level set of Hr is such that it is compact without singular points and contains 
both points (0(O),p0(O)) and (tt — (j){0),—p^{0)). In this case, the trajectory 
starting from (0(O),pJ(O)) with pj(0) > is periodic of period T and has a 
second crossing at the antipodal point tt — 0(0) after a time T/2. It is chased 
by a trajectory starting from {(j){0) , p^ (0)) where and p^ are roots of (fT5|) 
for (j> — (j){Q) with a time delay of T and reaches the antipodal point tt — (/)(0) 
at time T/2. 

Moreover, since the equations describing the evolutions of the remaining 
variables (r, 9) are invariant for the central symmetry: {(f>,p^) i-^ (tt — 0, — p^), 
we deduce that after half a period T/2, we have, for the two extremals starting 
respectively from (0(O),pJ(O)) and (0(O),p|^(O)) while (r(O),0(O) are identical, 
the relations: 

r+{T/2) = r-{T/2), 9+ {T/2) = 0- {T/2). (20) 

This generalizes the analysis of section 14.21 for Riemannian metrics on the two- 
sphere, except that we replace the reflectional symmetry by a central symmetry, 
leading to the fact that half a period instead of a quarter of period is necessary 
to construct the extremals. 
We proved: 

Proposition 16. If for fixed {pr^pe), the level set Hr = e is compact without 
singular point and has a central symmetry with respect to {4> — n/2,p^ = 0), then 
it contains a periodic trajectory of period T and if p^{Q) are distinct then 

we have two distinct extremal curves q^{t), <l~{t) starting from the same point 
and intersecting with the same length T/2 at a point such that <j){T /2) = tt— 0(0). 

Non-compact case: 

Non-compact level sets occur when p^ oo. Using (|13p . the relation (jj — 
gives the set S of solutions (j)s of 

(7.-r)^^-l = (21) 

and we must have |r — 7+I > 2. We deduce the following proposition: 

Proposition 17. //|r — 7+1 > 2 then we have extremal trajectories such that (f) 
is not periodic, i.e., — > 0, (j) (j)s and p^ — > ±00 when t —f +00 while — s- 
outside the equator. 
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Geometric interpretation: 

The two types of extremals are related to the Zermelo navigation problem 
of section 1331 

• If I7+ — r| < 2 then the system restricted to the two-sphere defines a 
Finsler geometric problem for which the extremals of the Grusin case are 
deformed into extremals described in proposition 1161 

• If I7+ — r| > 2 then we have two types of extremals: periodic extremals 
occur in a band near the equator and non periodic extremals occur when 
crossing the band around = 7r/4 where the current is maximal and 
asymptotic properties of proposition [17] correspond to the barrier phe- 
nomenon. 

This will be clarify by the pictures of the final section. 
Small time optimal analysis at the equator: 

We make a singularity analysis of the extremals near a point of the equator 
which can be identified to = 7r/2 and 6 — 0. Following the techniques of J], 
the main point is to construct a normal form in order to compute the small time 
optimal synthesis. 

Proposition 18. Near the equator, the small time optimal synthesis is given 
by the one near (x, y, z) — (0, 0, 0) for the system: 

(7+-r) 2 

X ^ 1 + y 

y = ir-j+)y + U2 (22) 
z = yui. 

Proof. The computation is straightforward. We set ■0 = n/2~(f> and near 7/1 = 0, 
we have 

f ~ -r + (r - ^+)ip^, - (r - 7+)^ - ""2, ^ ^ -ipui. 

This gives the normal form by setting x ~ —r/T, y ~ ip, z = 9 and us- 
ing a feedback transformation preserving \u\ < 1 and changing (mi,U2) into 

(-Ul,-U2). □ 

This normal form preserves the integrability property and the optimality 
analysis amounts to approximate the Hamiltonian Hr near an equatorial point. 
Moreover, by homogeneity, it describes the extremal behaviors not only near 
but also along an equatorial line identified to {t, 0, 0). 
From this normal form, we can deduce: 

Proposition 19. For the system i2S^) . near 0, every small time optimal trajec- 
tory is either 

1. A trajectory in the meridian plane z — Q with ui = such that: 

(i) //7+ — F < 0, an arc of the form BSB, i.e., a concatenation of a bang 
arc M2 = ±1, a singular arc and a bang arc, 

(ii) If ^+ — T > Q, a bang-bang arc BB, 

2. or a smooth extremal of order with pz ^ 0. 
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Proof. First of all, we analyze extremal curves contained in the plane z — 
for which ui = pz = 0. They correspond to the time-optimal control problem 
in the 2D-plane {x,y) with \u2\ < 1. They are already analyzed in [TJ using 
techniques of [13] , but the computed normal form reveals also the properties of 
the control. The line S : 1 1-^ {t,0) is a singular trajectory and the corresponding 
control is U2 = 0. Moreover, this line is time minimizing if 7+ — F < 0, while it 
is time maximizing if 7+ — F > 0. Other extremals are defined by U2 = sign(pj,). 
In case (i), a bang arc can have a connection with the singular arc, with a 
contact of order two with the switching surface. In contrast, the connection is 
not possible in case (ii). To complete the analysis, we observe that in case (i), 
every extremal curve which is of the form BSB is optimal, while in the case 
(ii), every extremal curve is bang-bang but for optimality, we have at most one 
switching. 

In order to conclude the analysis, we must prove that an extremal of order 
with pz cannot reach the singularity y = py = 0. The Hamiltonian Hr 
takes the form 

Hr - P. [1 + ^^^^^\ '] + Pyir - 7+)y + ^Pl+vlv^- 

Hr = e, with e = 0, 1, gives the condition: 

Px ^^j, ^ y^ +Py{T-j+)y+ ^pI+p'^zV'^ = 0. 
This is clearly not possible, taking the Taylor expansions of y, Py: 
y{t) = at + o{t), py{t) =ht + o{t) 

of an extremal of order 0, with pz ^ 0, reaching or departing from the singularity. 
The result is proved. □ 

As a corollary, we have the following result to determine the cut-locus in the 
case I7+ — F| < 2 where the time- minimal control problem is related to a Finsler 
problem. 

Corollary 1. Let q'^{t) and {t) be two distinct extremals of order with non- 
zero pe and intersecting with the same time T . Then they cannot he optimal 
beyond the intersecting point. 

Proof. The proof is standard. Assuming optimality beyond the intersecting 
point, we can construct a broken minimizer which is an extremal of order 
with non-zero pg. This contradicts proposition [TH □ 

5.2.3 Numerical computations 

We complete the analysis of section [5] by a series of numerical computations 
on extremal trajectories and on conjugate loci (see section [2]). All the compu- 
tations are done in the case of the time-minimum control problem but could 
be equivalently done for the energy minimization problem. We use numerical 
methods described in [7]: an extremal is computed and we calculate along this 
trajectory the conjugate points. 
Extremal trajectories: 
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We compute extremal trajectories for different values of the dissipative pa- 
rameters r and 7+. We represent in figures [T] and [2] in the coordinates (6*, </>) 
the projection of the trajectories on the sphere of radius 1. The variation of 
the radial coordinate p can be deduced straightforwardly and depends on the 
value of r and 7+. Extremals can also be plotted in the plane by using 

equation (fT^ . In figure [1] two trajectories intersecting with the same cost on 
the antipodal parallel are plotted. All other extremal trajectories in the case 
|r — 7+1 < 2 have the same qualitative behavior. For the Grusin model on the 
sphere, we also recall that the two trajectories which intersect on the cut locus 
have opposite initial values of p^, e and pe being fixed. This is no more the case 
when r ^ 7_|_ since these initial values depend now on T and 7+. The extremal 
trajectories for |r — 7+I > 2 are displayed in figure [3 We observe two types of 
trajectories, the periodic and the aperiodic ones. The periodic extremals have 
the same behavior as the ones of the case |r — 7+I < 2. When t — + -|-cxd, the 
aperiodic extremals have an asymptotic limit ips in (j) solution of the equation 
([2T|) . For given values of pe and e, this equation has two solutions symmetric 
with respect to the equator. As can be seen in figure [H we have found two 
types of aperiodic extremals: trajectories monotone in cf) which do not cross the 
parallel 4> — (t>s and trajectories passing by a maximum or a minimum in cf) 
different from 0s and crossing this parallel. These different types of behaviors 
can be determined by using equations (|17p and ([T8|l . 
Conjugate locus: 
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Figure 1: Extremal trajectories for F = 2.5 and 7+ = 2. Other parameters are 
taken to be ^^(0) = —I and 2.33, 0(0) ~ 7r/4, Pr = I and pg = 2. Dashed lines 
represent the equator and the antipodal parallel located at = 37r/4. 

We next determine numerically the conjugate locus for different set of pa- 
rameters. We restrict the discussion to the case jF — 7+I < 2. A similar study 
can be done for periodic trajectories if |F — 7+I > 2. For aperiodic trajectories, 
it can be shown that they are locally optimal in the sense that they do not 
have conjugate points for t G [0,-|-cxd[. Different conjugate loci are represented 
in figures O HI and \E\ both for the Grusin model and for a deformation of this 
model with the constraint jF — 7+I < 2. The conjugate locus of the Grusin 
model on the sphere is given in [6]. Here, in the case F = 7+, the drift vector 
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Figure 2: Extremal trajectories for T — 4.5 and 7+ = 2. Dashed lines represent 
the equator and the locus of the fixed points of the dynamics given by equation 
(PT|) . The solid line corresponds to the antipodal parallel. Numerical values of 
the parameters are taken to be (j){0) = 27r/5, pe = 8 and Pr{0) = 0.25. The 
different initial values of are -50, -10, 0, 2.637, 3, 5, 10 and 50. 

Fq being purely radial, the projection of the conjugate locus is the same as the 
conjugate locus of the Grusin model on the sphere. In particular, this projection 
is independent of the value of pr- Figure [3] displays both this conjugate locus 
and some extremal trajectories. Starting from the Grusin model, we can then 
modify the difference F — 7+ and determine the corresponding deformation of 
the conjugate locus. The evolution of the radial component being more com- 
plicated (the radial coordinate is no more decoupled from other coordinates if 
F 7^ 7+), the projection of the conjugate locus depends on the value oi Pr- A 
first comparison between the two cases is given by figure 2] where it can be seen 
that the global structure of the extremals is nearly the same. Figure [5] displays 
the projection of a conjugate locus on the plane {9, (p) corresponding to a par- 
ticular value oipr- The conjugate locus of the Grusin model has been added for 
the sake of comparison. We note that this locus is only slightly modified when 
F and vary. The trajectories of Fig. 0] are plotted in Fig. [5] up to the first 
conjugate point. 
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Figure 5: Projection of conjugate locus in solid line for pr — 0.5. The conjugate 
locus of the Grusin model corresponding to 7+ = F = 2 is represented in dashed 
lines. The horizontal dashed line indicates the position of the cut locus for the 
Grusin model. Dissipative parameters are taken to be F = 2.5 and 7+ = 2, pg 
is equal to 2. 
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Figure 6: Same as Fig. [3]but the trajectories are plotted up to the first conjugate 
point. 



